Exascale Ready Work-Optimal Matrix Inversion

In this thesis I present a new algorithm OPT for matrix inversion that builds on a matrix multiplication subroutine. It is combined of Strassen’s matrix inversion algorithm and Newton approximation. OPT overcomes the linear lower bound in parallel runtime of Strassen’s inversion algorithm and traditional Gaussian elimination without the log-factor more work of Newton approximation. In particular I prove, that given a work-optimal multiplication subroutine that runs in polylog time, OPT not only runs in polylog time, too, but furthermore only needs a constant factor more work than any work-optimal inversion algorithm. Additionally, I present a new stability result for Strassen’s matrix inversion algorithm combined with Newton approximation. As part of this thesis, I implemented OPT, Strassen’s inversion algorithm, and Newton approximation in a flexible test program along with a matrix container class optimized for this purpose. I describe the design of this implementation and the use and difficultys of BLAS for the matrix multiplication subroutine. In the experimental part, I compare the runtime of OPT to the routine included in the Intel Math Kernel Library (MKL) and observe its numerical stability. The constant factors on the amount of work of OPT show to be no larger than twice those of the MKL routine. As predicted, OPT shows to be very scaleable. Even on a computer with only eight cores it is already significantly faster than the MKL routine. Concerning numerical stability, OPT and Strassen’s algorithm do not live up to its bad reputation. Instead they produce results comparable to the MKL routine. I discover an unexpected instability of Newton approximation that makes it produce worse results than any other algorithm in the implementation. About this instability I present some further experiments.